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Summary 

The computer codes developed in this work provide a self-consistent set of thermo- 
dynamic and transport properties for equilibrium air for temperatures from 500 to 30000 K 
over a pressure range of 10 -4 to 10 2 atm. The properties computed include enthalpy, total 
specific heat, compressibility factor, viscosity, and the total values of thermal conductivity and 
Prandtl number. These properties are calculated through the use of temperature dependent 
curve-fits for seven constant values of pressure. A simple method of interpolation is demon- 
strated for intermediate values of pressure although other methods may be employed. The 
curve fits are based on mixture values calculated from an 11 -species air model. Individual 
species properties used in the mixture relations are obtained from a recent study by the 
present authors. A review and discussion of the sources and accuracy of the curve-fitted data 
used herein are given in NASA RP 1260. 

The computer subroutines developed to evaluate the various properties and described in 
this report are available from the NASA Computer Software Management and Information 
Center (COSMIC). 
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curve-fit coefficients for compressibility factor, Z, (eq. (10)) 
curve-fit coefficients for viscosity, p, (eq. (12)) 
curve-fit coefficients for total thermal conductivity, K 
curve-fit coefficients for total thermal conductivity, Pr , (eq. (17)) 
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total thermal conductivity, cal /cm -sec -K 

total Prandtl number, C p \i/K 

pressure, atm 

temperature, K 

compressibility factor 

viscosity, poise 

= (771000) in Eqs. (11), (13), and (17) 

= In (7710,000) in Eqs. (6), (9), and (16) 


Analysis 

It is desirable to have curve fits to the various equilibrium air properties as opposed to 
tabulated values since curve fits generally permit more efficient computation for flow-field 
analyses. In addition, for accurate calculations, it is preferable that these curve fits be obtained 
to values which are computed in a self-consistent manner from the same set of data and other 
physical properties. The present work follows these guidelines closely. In this report, details 
of the curve-fits and associated computer routines developed to evaluate the equilibrium pro- 
perties are presented along with a brief physical overview. The reader should consult the 
companion report (Ref. 1) for a more detailed discussion of the sources and accuracy of the 
methods employed to calculate the original equilibrium air data. The computer codes 
developed in this work to evaluate the thermodynamic and transport properties of equilibrium 
air may be obtained for a fee from: 

COSMIC 

Computer Services Annex 

University of Georgia 

Athens, GA 30602 

(404) 542-3265 

Request the codes by the designation LAR 14760. 

These codes were written in the form of FORTRAN subroutines for easy adaption to 
existing programs. The subroutines are well commented and can be easily modified to suit 
the user’s needs. 

The present equilibrium air properties are based on calculations using an 11-species ( N 2 , 
O 2 > O , NO, NO + , e , N + , 0 + , N++, and 0 ++ ) air model and which includes the single 
and double ionization effects encountered at high temperatures. The curve fits provided here 
can only be used for equilibrium air calculations without ablation or injection of foreign pro- 
ducts, i.e., with constant elemental composition through the flowfield. Values of 0.24 and 
0.76 are employed in this work for the elemental mass composition of oxygen and nitrogen, 
respectively. In addition, these properties assume charge neutrality for ionized air and, there- 
fore, neglect the associated electric field effects. This assumption can generally be made for 
an ionized gas mixture in thermochemical equilibrium. 

The present equilibrium thermodynamic and transport properties are for temperatures 
from 500 to 30000 K over a pressure range of 10 -4 to 10 2 atm. Included in these properties 
are enthalpy, total specific heat, compressibility factor, viscosity, and the total values of ther- 
mal conductivity and Prandtl number. Curve fits of these properties are based on the use of 
temperature dependent polynomial curve-fits for discrete values of pressure. Seven pressure 
intervals from 1(T 4 to 10 2 atm in constant steps by factors of ten were considered. A simple 
linear interpolation in the logarithm of the pressure is employed for intermediate values. 
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Other interpolation schemes may be developed to cover these intermediate areas. 

The order of the curve-fits polynomials and the number of temperature intervals for the 
curve-fits have been selected so as to achieve an overall 5 percent accuracy of the curve-fit 
values in relation to the computed values along the constant pressure lines. It was also desir- 
able to employ the minimum number of temperature intervals (or ranges) and the lowest order 
polynomials possible along the constant pressure lines while achieving this accuracy. In the 
final result, the entire temperature range of 500 < T < 30000 K is broken down into two or 
more intervals using these criteria. To ensure a smooth variation of the curve-fitted properties 
over the complete temperature range, values of the curve-fit coefficients are linearly averaged 
at the temperature interval boundaries over a specified overlap region. The overlap range is 
equally spread on the two sides at the temperature boundaries and is 500 K for temperatures 
less than 20000 K and 1000 K for higher temperatures. The averaging procedure for all 
coefficients in the overlap region can be written as 


A* = ( \-a)A\+aAV 

where 


( 1 ) 


T -Tk 



( 2 ) 


Here, and Tjj are the lower and upper temperature limits, respectively, of the overlap 
region. The quantity is one of the curve-fit coefficients for flowfield property "<J>". The 
superscripts U and L on denote coefficients of the curve fits from the upper and lower 
sides of the temperature boundaries, respectively. With these definitions, all curve-fit 
coefficients A ^ are replaced by from Eq. (1) when the temperature (at which the property 
is to be evaluated) is within an overlap region. 


Since the curve fits are obtained for a range of constant values of pressure (KT^p^lO 2 ) 
in steps of a factor of 10, a logarithmic interpolation procedure was used to obtain values of 
the flowfield properties at intermediate pressures. That is, a property "<J>* " is interpolated at 
pressure p (for a specified temperature) from the following relation: 


log<(>* 


log ^2 ~ iogfrl) 

(logp 2 - logPi, 


(logp* - logp t ) + log <(>j 


(3) 


where <j> 2 and 4>j are the values of property <|) evaluated for a specified temperature T at pres- 
sures p 2 and p j, respectively, such that p { <p* <,p 2 . 

Using these relations and the curve-fits described subsequently, the interpolated values of 
viscosity (p), specific heat at constant pressure ( C p ), thermal conductivity ( K ), and Prandtl 
number ( Pr ) are within 5 percent of the computed values for temperatures less than 2500 K 
over the entire pressure range. The error in the interpolated values of p, C p , K, and Pr 
increases at higher temperatures. A maximum error of about 20 percent is obtained in the 
interpolated values of C p and Pr at temperatures between 5000 K and 8000 K, whereas this 
error is about 10 percent or less for K for temperatures between 2500 K and 8000 K. For 
temperatures greater than 8000 K, the maximum error in the interpolated values is believed to 
be less than 20 percent for C p ,K , and Pr, and much less for p. Such errors in the interpo- 
lated values of C , K, and Pr are compounded by the presence of the reaction component in 
these properties. The interpolation errors verified in the present study and expected using the 
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logarithmic scheme are summarized in the following table. 

Interpolation Errors (for P=0.03, 0.3, and 3.0 atm) 


Temperature Range, K 

Property 

Error 

500 < T < 15000 

V 

< 8% 

500 < T < 15000 

h 

< 5% 

500 < T < 15000 

Z 

< 5% 

500 < T < 2500 

C p 

< 5% 

2500 <T <, 5000 

4 

< 10% 

5000 < T < 15000 

Cp 

< 20% 

500 < T < 2500 

K 

< 5% 

2500 < T < 15000 

K 

< 10% 

500 < T < 2500 

Pr 

< 5% 

2500 < T < 5000 

Pr 

< 10% 

5000 <T< 15000 

Pr 

< 20% 


Finer pressure intervals (than those used here) may improve the interpolation accuracies by 
capturing the initiation of various reactions better. Additionally, a different approach to inter- 
polating for properties at the intermediate pressures may be developed which will improve the 
accuracy. However, even with the interpolation scheme and pressure intervals employed 
presently, flowfield quantities of interest such as the skin-friction drag, surface heat transfer, 
and the aerodynamic coefficients can be computed for cooled surfaces with fair degree of 
accuracy. 

The curve fits provided here are in terms of temperature. These may easily be converted 
to those in terms of enthalpy by using the relation between temperature and enthalpy provided 
by the curve fits given in the next section. 

Enthalpy 

The specific enthalpy for equilibrium air can be computed from following curve-fit 
expression to obtain enthalpy as a function of temperature at constant pressure 



h =exp (A h x 4 +B h x 3 + C h x 2 +D h x+E h ) 

(4) 

or. 


In h = A h % 4 +B h x 3 +C h x 2 +D h x + E h 

(5) 

where 


X = In (7710,000) 

(6) 


The coefficients Ay, , B y, , C y, , Dy, , and have been evaluated from the computed values 
by using a least squares curve-fit technique and are given in Table I. These coefficients are 
given for several temperature ranges at a specified pressure. As stated previously, the number 
of temperature ranges is established so as to obtain a smooth and accurate curve fit to the 
computed values and generally, a larger number of ranges is needed at lower pressures. The 
accuracy of the curve fit expression (4) is good to within 4 percent of the computed values 
(Ref. 1) for the entire temperature range (500 K <, T <, 30000 K) range along the constant 
pressure lines. 
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Total Specific Heat At Constant Pressure 

The total specific heat at constant pressure can be computed from the following curve-fit 
expression 


C p = exp(A Cp x 4 + B Cp x 3 + C Cp x 2 +D Cp X + E Cp ) 

(7) 

or. 



where 

C p = A c p X 4 + B c p X 3 + C Cp X 2 + D Cp X + E Cp 

(8) 


X = In (7710,000) 

(9) 


The tabulated values of coefficients A Cp , B c p , C Cp , D Cp , E Cf are listed in Table II. The accu- 
racy of the curve-fit values is within 4 percent of the computed values (Ref. 1) along constant 
pressure lines for the entire temperature range considered here. 


Compressibility Factor 

The compressibility factor, Z, (which is equal to the ratio of the molecular weight of 
undissociated air to the mean molecular weight) can be obtained by employing the following 
curve-fit expression at constant pressure: 

Z = A z +B z x + C z x 2 +D z x 2 + E z x a (10) 

where 


T 

1000 


( 11 ) 


The coefficients appearing in Eq. (10) are tabulated in Table III. The curve fit values are 
within 3 percent of the computed values (Ref. 1) for values along the constant pressure lines. 


Viscosity 

The following expression has been used to curve-fit the viscosity values computed in 
Ref. 1, 


|i = A il + B il x + C il X 2 +D il x 3 + E v X A+F i iX 5 (12) 

where 


T 

1000 


(13) 


The polynomial coefficients of Eq. (12) are given in Table IV. In general, the accuracy of the 
curve fits along constant pressure lines is good to within 4 percent of the value computed in 
Ref. 1. The accuracy improves with increasing pressure, similar to the other thermodynamic 
and transport properties. 


Total Thermal Conductivity 

The total thermal conductivity, K , for equilibrium air is obtained from the curve-fit 
expression 

K = cxp(A K x 4 + B K x i + C K x 2 + D K x + E K ) (14) 

or 
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In K - A k x 4 + B k x 3 + C k x 2 + D k x + E k (15) 

where 

X = In (7710,000) (16) 

The coefficients appearing in Eqs. (14) and (15) are given in Table V. The accuracy of the 
curve fit is within 4 percent of the computed values along the constant pressure lines for the 
entire temperature ranges. 


Total Prandtl Number 

The computed values of the total Prandtl number have been curve fit by employing the 
relation 

Pr = A Pr +Bp r x + Cp r x 2 + D Pr x 3 + E Pr x 4 +Fp r X 5 (17) 

with 


T 

1000 


(18) 


The polynomial coefficients of Eq. (17) are provided in Table VI. The accuracy of the curve- 
fits is within 4 percent of the values computed in Ref. 1 for the constant pressure and variable 
temperature ranges considered here. 


User Instructions 

The computer codes developed in this work to evaluate the thermodynamic and transport 
properties of equilibrium air were written in the form of FORTRAN subroutines for easy 
adaption to existing programs. The subroutines are commented and can be easily modified to 
suit the user’s needs. In an attempt to maintain generality, a total of six separate subroutines 
are available for use. These six routines are used for the evaluation of the following proper- 
ties: 

1) specific enthalpy, h, kcal/g 

2) total specific heat at constant pressure, C p , cal/g-K 

3) compressibility factor, Z 

4) viscosity, 1 1 , poise 

5) total thermal conductivity, K , call cm -sec -K 

6) total Prandtl number, Pr 

Each routine evaluates one of the curve-fits described above at a given (input) temperature 
and pressure and returns a dimensional quantity for the associated property. It should be 
noted that all of the curve-fit coefficients used in the present work are listed in the form of 
DATA statements in the subroutines where they are used. In this form, it would be a simple 
task to incorporate the coefficient data into computer codes developed independently of those 
described in this report. The routines could also be modified to return the properties for 
arrays of temperatures and pressures for improved efficiency in some applications. 

Subroutine Description and Use 

The calling syntax of the six subroutines developed in the present work are outlined in 
the following sections. Each subroutine was named to reflect the property computed as shown 
in the following list: 
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1) ENTHLPY - specific enthalpy 

2) SPECIFC - total specific heat at constant pressure 

3) COMPRES - compressibility factor 

4) VISCSTY - viscosity 

5) CONDUCT - total thermal conductivity 

6) PRANDTL - total Prandtl number 

In the following descriptions, parameter lists are described for each subroutine. 

Specific Enthalpy 

- ENTHLPY 

Purpose: Compute specific enthalpy at a given temperature and pressure 

Use: CALL ENTHLPY(P,T,H) 

Input: P = Pressure, atm (1CT 4 < P < 10 2 ) 

T = Temperature, K (500 < T < 30000) 

Output: H = h , specific enthalpy, kcal /g 

Specific Heat 

- SPECIFC 

Purpose: Compute the specific heat at constant pressure at a given temperature and 

pressure 

Use: CALL SPECIFC(P,T,CP) 

Input: P = Pressure, atm (10 -4 < P <, 10 2 ) 

T = Temperature, K (500 < T < 30000) 

Output: CP = C p , specific heat at constant pressure, g/cm-s 

Compressibility Factor 

- COMPRES 

Purpose: Compute the compressibility factor at a given temperature and pressure 

Use: CALL COMPRES(P,T,Z) 

Input: P = Pressure, atm (10 -4 < P £ 10 2 ) 

T = Temperature, K (500 < T < 30000) 

Z = Z, compressibility factor 


Output: 
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Viscosity 


- VISCSTY 

Purpose: Compute the viscosity at a given temperature and pressure 

Use: CALL VISCSTY(P,T,AMU) 

Input: P = Pressure, atm (lO^ 4 < P <, 10 2 ) 

T = Temperature, K (500 < T < 30000) 

Output: AMU = p, viscosity, poise 

Thermal Conductivity 


- CONDUCT 

Purpose: Compute the total thermal conductivity at a given temperature and pressure 

Use: CALL CONDUCT(P,T,AK) 

Input: P = Pressure, atm (10" 4 < P < 10 2 ) 

T = Temperature, K (500 <, T < 30000) 

Output: AK = K , total thermal conductivity, cal /cm— sec -K 

Prandtl Number 


- PRANDTL 

Purpose: Compute the total Prandtl number at a given temperature and pressure 

Use: CALL PRANDTL(P,T,PR) 

Input: P = Pressure, atm (10 -4 £ P < 10 2 ) 

T = Temperature, K (500 < T < 30000) 

Output: PR = Pr , total Prandtl number 


Environmental Characteristics 

The computer codes described in this report were developed on a SUN SPARCstation 2 
computer running the SunOS 4.0 operating system. The codes are written in standard FOR- 
TRAN 77 language and should run without modification on other machines. All routines are 
written with double precision (64 bit) reals for accuracy in evaluating the relatively high order 
polynomials employed. The memory sizes of the various routines are tabulated in the follow- 
ing table: 
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Memory Size of Subroutines 


Routine 

Size, bytes 

ENTHLPY 

3853 

SPECIFC 

1395 

COMPRES 

1543 

VISCSTY 

3023 

CONDUCT 

3039 

PRANDTL 

3039 


The codes were timed on a Sun SPARCstation 2 with global optimization invoked for 
the FORTRAN compiler (i.e., f77 -03). Times were monitored for 10 4 calls to each subrou- 
tine for an intermediate pressure (3 atm ), and the average time for a single call was computed 
as shown in the following table: 

CPU Time for Single Subroutine Calls 


Routine 

Time, s 

ENTHLPY 

8.7 x 10~ b 

SPECIFC 

9.5 x 10" 5 

COMPRES 

7.2 x 10~ 5 

VISCSTY 

7.0 x 10- 5 

CONDUCT 

9.2 x 1(T 5 

PRANDTL 

7.5 x 10" 5 


Implementation Instructions 

The compilation command for the FORTRAN compiler on the Sun SPARCstation 2 sys- 
tem for the subroutines described in this guide is 

P7 -03 -c filename f 

where filename / is the specific program name from the list of names described earlier (i.e., 
ENTHLPY.SPECIFC,. ...PRANDTL). The object file from the compilation step above would 
be linked to a user program by the following command: 

P7 -ffpa -03 -o exec userprogram.o filename.o 

where exec is the executable file to be created, merprogram.o is the object file of the user’s 
program after compilation, and filename.o is the object file from the subroutine compilation. 
In the compilation and linkage steps above, the command line parameter -03 is used to 
invoke global optimization which may be omitted if desired. 

A driver program was written along with the subroutines described thus far which is 
used to exercise each of the subroutines for typical conditions. This driver program is useful 
for testing the subroutines and could be modified to generate a useful collection of tabulated 
properties. Output from the driver program is listed in the next section. 

The sequence of commands listed below is used to execute the driver program on the 
Sun SPARCstation 2 system. 
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J77 -03 -o DRIVER DRIVER f ENTHLPY.f SPEClFC.f COMPRES.f 
VISCSTY.f CONDUCT./ PR ANDTL.f 

DRIVER 

Sample Input/Output 

The driver routine described above is used to exercise each of the subroutines for typical 
input parameters. This section lists the output from the driver program which is generated 
directly using the output from the various subroutines. The driver routine calls all six subrou- 
tines to evaluate the equilibrium air properties for a range of temperatures (T) from 1000 to 
25000 K for a constant pressure of 1 atm . The output from each subroutine is gathered by 
the driver program for all temperatures and is then printed in tabular form. The following 
table represents the output from the driver program. 

Output from the DRIVER program 


T, K 

H, KCAL/G 

CP, CAL/G-K 

Z 

MU, POISE 

K, CAL/CM-S-K 

PR 

1000.0 

0.25541 

0.27521 

0.99354 

0 . 4 1558E-03 

0.1634 0E-O3 

0.70378 

2000.0 

0.54646 

0.31575 

0.99734 

0 . 65842E-03 

0 . 3085 9E-03 

0.68567 

3000.0 

0.99546 

0.63713 

1.0438 

0 . 85660E-03 

0 . 98728E-03 

0.58220 

4000,0 

1.9005 

0.77051 

1.1258 

0. 10596E-02 

0. 14361E-02 

0.59570 

5000.0 

2.4441 

0.66248 

1.2233 

0 . 12962E-02 

0. 14556E-02 

0.59818 

6000.0 

3.7070 

1.8092 

1.3027 

0.15742E-02 

0 . 4 8761E-02 

0.56958 

7000.0 

6.2717 

3.2664 

1.6070 

0.18809E-02 

0 . 82702E-02 

0.73562 

8000.0 

9.2197 

1.8646 

*r 

CT\ 

r- 

CD 

0 . 2 1 824E-02 

0 . 4 64 97E-02 

0.90105 

9000.0 

10.532 

1.0399 

1.9864 

0 . 2424 IE-02 

0 . 2833 6E-02 

0.88022 

10000.0 

11.497 

1.0966 

2.0334 

0 . 25305E-02 

0 . 2 9954E-02 

0 . 93840 

11000.0 

12.859 

1.7093 

2.1009 

0 . 24050E-02 

0.37364E-02 

1.0738 

12000.0 

15.018 

2.5369 

2 .2240 

0 . 20425E-02 

0 . 4 662 8E-02 

1.1410 

13000.0 

18.161 

3.7196 

2.4229 

0 . 1 5550E-02 

0 .55023E-02 

CD 

in 

o 

14000.0 

22.310 

4.8176 

2.6949 

0 . 10822E-02 

0 . 6007 3E-02 

0.85147 

15000.0 

27.242 

5.1329 

3.0149 

0.70535E-03 

0.61244E-0 2 

0.58859 

16000.0 

32.374 

4.5418 

3.3395 

0 . 4 4 773E-03 

0 .59685E-02 

0.34731 

17000.0 

36.752 

3.5076 

3.5748 

0.29438E-03 

0.57214E-02 

0.18211 

18000.0 

39.583 

2.5291 

3.735B 

0.2U55E-03 

0 . 5554 9E-02 

0 . 94 5 96E-03 

19000.0 

41.775 

1.8305 

3.8390 

0. 16621E-03 

0 . 561 52E-02 

0 .53406E-01 

20000.0 

43.262 

1.4279 

3.9002 

0 . 14582E-03 

0 . 58 621E-02 

0.35340E-01 

21000.0 

44.453 

1.2553 

3.9350 

0 . 13808E-03 

0 . 623 91E-02 

0.27895E-01 

22000.0 

45.666 

1.3206 

3.9587 

0 . 13708E-03 

0.67285E-02 

0.261 12E-0 1 

23000.0 

47.140 

1.5895 

3.9857 

0 . 13809E-03 

0 . 72978E-02 

0 . 28862E-0 1 

24000.0 

49.060 

2.0770 

4.0277 

0 . 13776E-03 

0.79202E-02 

0 .35793E-01 

25000.0 

51.579 

2.8165 

4.0899 

0 . 13405E-03 

0.85734E-02 

0.441 90E-0 1 
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Table I. Curve-fit coefficients for enthalpy* of equilibrium air** (500 K < T < 30000 K)t 
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Table II. Curve-fit coefficients for specific heat* of equilibrium air (500 K < T < 30000 I\) 
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Table Ml. Curve-fit coefficients for compressibility of equilibrium air (500 K < T < 30000 K) 
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Table IV. Curve-fit coefficients for viscosity* of equilibrium air (500 K < T < 30000 K) 
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Tcihle VI. Curve- fit coefficients for Prandtl number of equilibrium air (500 K < T < 30000 I\) 
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